Dynamical identification of open quantum systems

I propose a quantum-trajectories approach to parametric identification of the effective Hamiltonian for a Markovian open quantum system, and discuss an application motivated by recent experiments in cavity quantum electrodynamics. This example illustrates a strategy for quantum parameter estimation that efficiently utilizes the information carried by correlations between measurements distributed in time. The rapidly developing theory of quantum parameter estimation (QPE) focuses on the design of optimal measurement strategies for extracting information about c-number parameters E θ that characterize a given quantum system. While much progress has been made in applying QPE to the parametric identification of prepared quantum states, little or no attention has yet been paid to the problem of estimating parameters that characterize a dynamical quantum system. In this letter I consider the latter aspect of QPE within the context of quantum optics, and describe a quantum-trajectories method for estimating parameters that appear in the effective Hamiltonian for a Markovian open quantum system. Recent theoretical work in QPE [1–3] has focused on a paradigm in which an experimenter is provided with one or more copies of a quantum state ρ0 drawn from a single-parameter family ρ(θ), and is asked to determine the value θ0 such that ρ0 = ρ(θ0). The experimenter knows the form of ρ(θ) and can make arbitrary measurements on the states s/he is given, but does not know the value of θ0 a priori. In such situations one can actually derive a mathematical representation of the optimal quantum measurement for the purpose of estimating θ0, and optimize over all possible statistical reductions of the measurement results as well [4]. Accordingly, there exists a generalized version of the Cramér–Rao inequality [5] that establishes a fundamental bound on the rate of convergence for estimators based on repeated measurements whose marginal statistics are fully determined by a single, unchanging density matrix [6]. My purpose here is to consider a related but distinct aspect of QPE, namely the estimation of parameters appearing in the equations of motion that govern the time-evolution of a quantum system. In this paradigm, which is closely related to that of classical system identification [7], the hypothetical experimenter wishes to determine which system model H0 in a parametrized family H(E θ) best accounts for the dynamical behaviour of a given quantum system. In contrast to the conventional QPE scenario described above, the statistics of multiple successive measurements made on a dynamical system cannot necessarily be derived from any single density matrix. Roughly speaking, this is because every measurement disturbs the state of the system [8] in a manner that depends on the † E-mail address: [email protected] 1355-5111/96/061103+06$19.50 c © 1996 IOP Publishing Ltd 1103 1104 Letter to the Editor measurement outcome, and because the evolution of the system state between measurements depends on the full details of the system’s equations of motion. The effects of repeated measurements on otherwise unitary quantum evolution have been extensively studied in quantum optics, with regard to the dynamics of open quantum systems [9, 10]. The configuration most often treated by such work is that of a small, ‘encapsulated’ quantum system having one or more well defined input/output channels associated with its coupling to the physical environment. This picture naturally suggests a paradigm in which an experimenter attempts to parametrically identify the Hamiltonian of the encapsulated quantum system by examining the response of the output channels to driving stimuli applied to the input channels†. The task of quantum system identification may then be equated with that of computing the relative likelihood of an observed sequence of measurement results 4∗ as a function of the parameter set E θ , given the external driving conditions imposed by the experimenter. To the extent that the environmental couplings for the system are known, quantum trajectory theory [10, 12] suggests a simple method for the computation, which I discuss below. Having a likelihood function f (E θ | 4∗), one can use maximum-likelihood or Bayesian principles [13, 7] to estimate the parameters E θ . Note that it should generally be possible for the experimenter to determine optimal driving conditions that make the system response maximally sensitive to the values of E θ , or indeed to adaptively change the driving conditions as the estimation starts to converge [3]. To illustrate quantum system identification in a concrete setting, let us focus on an example with relevance to current experiments in cavity quantum electrodynamics (QED)—a single two-level atom placed within the mode volume of a driven, highfinesse optical cavity [14, 10]. The strength of the coherent coupling between atom and cavity mode is parametrized by the vacuum Rabi frequency g, whose value depends on the spatial position of the atom within the cavity. For a Fabry–Perot resonator g(E r) = g0 cos(2πx/λ) exp [−(y2 + z2)/w2 0], where x is the coordinate along the cavity axis and w0 is the Gaussian waist of the TEM00 resonator mode. The specific task I shall consider is that of estimating g ∈ [0, g0], which I suppose to be unknown because the atomic position is not known. The measurement procedure will simply be to monitor the arrival-time statistics of photons emitted by the atom–cavity system for a fixed cavity driving field. For the purposes of this discussion I shall not explicitly treat the atom’s external degrees of freedom, imagining that they are fixed by an RF Paul trap or similar confining mechanism [15]. However, note that the correlation of g with the atomic position operator implies that ‘online’ estimation of g for an untrapped atom drifting through a cavity could be viewed as a time-distributed quantum measurement of the position of a free mass [16, 17]. For a gedanken-experiment in which the cavity is driven by a resonant CW probe laser and both the atomic fluoresence and cavity emission are continuously monitored by perfect photon-counting detectors, the evolution of the conditional state-vector between photodetection events satisfies the effective Schrödinger equation (h̄ = 1) |ψc(t + dt)〉 = e−iHeffdt |ψc(t)〉 (1) Heff(g) = ig ( aσ+ − a†σ− ) + i (a − a†) − iκa†a − iγ⊥σ+σ. (2) This interaction-picture expression for Heff(g) is valid under the rotating-wave and electricdipole approximations, and for identical atomic/cavity/probe-laser frequencies. Here κ is the field decay rate of the cavity, γ⊥ is the dipole decay rate of the atom, and represents the strength of the coherent driving field. The jump operator associated with the detection † Wiseman and Milburn [11] have derived closed-loop master equations to describe the effects of feeding an observed output channel back to a controlled system input. Letter to the Editor 1105 of photons spontaneously emitted by the atom is ĉ0 ≡ √ 2γ⊥σ−, and ĉ1 ≡ √ 2κa is the jump operator associated with the detection of photons leaking through the cavity mirrors. By registering the origins {j1, . . . , jn} (= 0 for spontaneous emission or 1 for cavity decays) and arrival times {t1, . . . , tn} of every photon emitted by the atom–cavity system in response to the CW driving field during an observation interval [t0, tf ], the hypothetical experimenter accumulates a classical record 4∗ ≡ (t0, tf , {ji, ti}) of the stochastic evolution of the system state. Assuming a uniform prior distribution on g, the likelihood function f (g | 4∗) then simply corresponds to a normalized version of the exclusive probability density [10, 12] p ( 4∗| g) = Tr  Ueff ( tf , tn ; g ) ĉjn Ueff (tn, tn−1 ; g) · · · × ĉj1Ueff (t1, t0 ; g) ρ(t0) U † eff (t1, t0 ; g) ĉ† j1 × · · · U † eff (tn, tn−1; g) ĉ† jn U † eff ( tf , tn; g )  (3) viewed as a function of g rather than 4∗. Accordingly, the maximum-likelihood estimate (MLE) of g is obtained by computing the value of g which maximizes (3) with 4∗ fixed by the observed data. Here Ueff ( t ′, t ; g) is the evolution operator from time t to t ′ associated with the effective Hamiltonian Heff(g) defined in equation (2). In order to numerically demonstrate quantum system identification using (3), I have generated a set of classical records by quantum Monte Carlo simluation [21] of a driven atom–cavity system with (g0, γ⊥, κ)/2π = (57, 2.5, 30) MHz, and three different powers for the driving field = {24, 34, 44.3}. While this value for g0 is larger than what has been achieved experimentally [22], it should certainly be within the reach of works in progress. For the simulations I chose an arbitrary atomic position such that g(E r) = 45 MHz, and generated classical records with an observation time of 1 μs each. Figure 1(a) illustrates the stochastic time-evolution of the mean intracavity photon number, taken from typical Monte Carlo data sets for each of the three values of . The photocount statistics are clearly super-Poissonian, and the simulated data show that quantum jumps often occur at a local rate that greatly exceeds the rate at which the system regresses to steady state. Figure 1. (a) Time-evolution of the mean intracavity photon number 〈a†a〉 in individual trajectories. (i) is for = 24, (ii) is for = 34, and (iii) is for = 44.3. (b) Corresponding stochastic evolution of the (normalized) likelihood function f (g|4∗) and corresponding MLE in one quantum trajectory with driving field amplitude = 34. The surface height indicates the relative probability of g ∈ [35, 57], with a resolution of 1. The ‘true’ value of g corresponds to 45. Note that the likelihood function is updated each time a photon is detected, so that the timelike coordinate in this surface plot corresponds to jump number rather than absolute time. 1106 Letter to the Editor Figure 2. (a) Standard deviation of the maximum-likelihood estimator for g as a function of absolute time (+: = 24, ◦: = 34, ×: = 44.3). (b) Histogram showing the evolution (in absolute time) of the sampling distribution for the MLE of g, representing 150 simulations with = 44.3. For each Monte Carlo trajectory, an identification routine based on (3) was used to compute the stochastic time-evolution of f (g | 4∗), as well as the corresponding MLE. Figure 1(b) shows one typical data set with = 34, starting from the inital estimate made after only one photodetection event and updated after each subsequent photodetection. Figure 2(a) indicates the ensemble-averaged convergence of the MLE for g, based on 2000 simulations for = 24, 300 simulations for = 34, and 150 simulations for = 44.3. A histogram representing the time-evolution of the MLE sampling distribution is given in figure 2(b), for the case of = 44.3. With this driving field, ∼ 1% accuracy in estimation of g is obtained in 1 μs observation time (∼ 600 quantum jumps). It is important to note that the QPE procedure described above automatically makes efficient use of any information about g that is contained in higher-order correlations of the classical record of counting times. Of course, not every open quantum system will generate significant correlations of this type. In the scenario discussed above, for example, the photon stream emitted by the atom–cavity system would become nearly Poissonian in the limit of either weak excitation (correlated pairs of photons become rare) or of weak coupling g κ, γ⊥ (correlations become weak). Two critical conditions for correlations to be strongly evident in individual classical records are that the mean time between counts must be comparable to or less than the system regression time [19], and that the system dynamics must be significantly altered by the loss of a single quantum of excitation [23]. For systems not satisfying these criteria, the methods described above offer no real advantage over statistical estimators based on only the steady-state density matrix obtained by solving the master equation associated with (2). Accordingly, optimal parameter estimation in such systems can be formulated within the paradigm of conventional QPE. Judging from the trend shown in figure 2(a), however, it certainly seems that the information rate on parameters in strongly coupled systems can be significantly larger in the strong-driving regime [24]† than in the weak-field regime. † Indeed, one would like to quantify the information rate associated with the observation of classical records for different driving conditions. This would entail computing the Fisher information on E θ contained in the stochastic ensemble of records 4. Note that for a fixed effective Hamiltonian, the Fisher information rate on a given parameter could be optimized by varying the basis in which the output channels are measured. Letter to the Editor 1107 In closing, let me note that a straightforward extension of the above method would allow the identification of non-stationary Hamiltonians in which the parameters E θ(t) vary slowly compared to the timescale for convergence of the corresponding statistical estimators. Relative to the example discussed above, a recent cavity–QED experiment incorporating a laser-cooled atomic source [25] has demonstrated the practical feasibility of achieving this separation of timescales. It seems reasonable to hope that the methods proposed above could be utilized in future experimental work to track variations in g associated with the motion of an individual atom through the mode-volume of a high-finesse optical cavity [18]. A digital signal processor implementing such a procedure could be used as the state observer in a ‘semiclassical’ feedback control loop designed to confine and cool the atom’s centre-of-mass motion. I wish to acknowledge invaluable conversations with H J Kimble, G J Milburn, S M Tan, Q A Turchette, D W Vernooy and P Zoller, as well as the support of a National Defense Science and Engineering Graduate Fellowship.